Unraveling Self-similar Behavior and Periodicity in the Bouncing Ball through Wavelets 
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The non- stationary dynamics of a bouncing ball comprising of both periodic as well as chaotic behavior, is 
studied through wavelet transform. The multiscale characterization of the time series discloses clear signature 
of self- similarity and periodicity The scale dependent variable window size of the wavelets aptly captures both 
the transients and non- stationary periodic behavior The optimal time-frequency localization of the continuous 
Morlet wavelet is found to delineate the scales corresponding to different periodic modulations. Self-similar 
behavior is quantified by the generalized Hurst exponent, obtained through both wavelet based multifractal 
detrended fluctuation analysis and Fourier methods. The discrete Daubechies basis set has been applied for 
detrending the temporal behavior to reveal the multifractal character underlying the dynamics. 

Keywords: Wavelet based multifractal detrended fluctuation analysis (WBMFDFA), Morlet wavelet, Daubechies4 (Db4), 
Scaling exponent, Fractal dimension, Hurst & Lyapunov exponents. 



J2 

i 

B 



CO 

o 
m 



X 



INTRODUCTION 

Driven non-linear systems are well known to exhibit com- 
plex dynamics, ranging from stable to unstable periodic mo- 
tions and chaotic behavior, under different conditions [1-4]. 
A number of simple systems e.g., coupled pendulum, electri- 
cal circuits can reveal extremely complex dynamics [1]. Well 
known systems such as bilinear oscillator [5], billiards [6-7], 
Chua circuit [8-9] have been studied quite extensively in the 
literature. More complex systems like impact oscillators re- 
veal highly complex dynamical behavior, as certain periodic 
orbit of the system goes through grazing and enters the im- 
pacting regime [5]. Quantization of the Sinai billiard, a two 
dimensional stochastic and ergodic system, has been studied 
extensively. Bifurcations in the periodic orbit of classical dy- 
namical systems e.g., quartic oscillator, have been shown to 
have effect on the quantal dynamics [10]. All these make ex- 
ploration of classical non-linear systems, an exciting area of 
ongoing research. 

The non-linearity in these systems can arise from non-linear 
interactions or from imposed boundary conditions. In linear 
cases, the driving force can induce modulations to commen- 
surate with the driving frequency, whereas in the case of non- 
linearity, other frequencies can manifest. Here, we explore 
the rich dynamics of the driven ' Bouncing ball ' system and 
characterize its complex temporal behavior through Fourier 
method and wavelet based approach. This system has been 
studied earlier due to its interesting dynamics. Tufillaro et al., 
illustrated the period doubling route to chaos in such systems 
in the presence of dissipation [11-12]. They used frequency 
and amplitude as controlling parameters to demonstrate the 
curves in parameter space that separate periodic orbits. Heck 
et al., used high speed camera for observing the dynamics 
which was then modeled numerically [13]. Kini et al., con- 
figured a computer controlled system in order to ascertain the 
state of the bouncing ball system more accurately [14]. A. 
Okninski et al. have studied numerically as well as analyt- 
ically the dynamics of a bouncing ball moving in a gravita- 



tional field and colliding with a vertically moving table with 
constant velocity [15]. 

We identify and locate the time varying periodic com- 
ponents, as well as quantify the chaotic dynamics through 
wavelets. This driven system shows non- stationary behavior 
and unstable periodic orbits. It also exhibits chaotic dynamics 
under appropriate parametric conditions. We make use of both 
continuous and discrete wavelets to characterize the periodic 
as well as self similar behavior of this system [16]. The pe- 
riodic components are ascertained by isolating the variations 
in both temporal and frequency domains. Morlet wavelet, a 
product of a variable Gaussian window and a sinusoidal func- 
tion, is used to identify local periodic modulations. Wavelet 
transform splits the input time series into components that de- 
pend on position and scale. One then characterizes the varia- 
tions in it, by changing the scale for a particular location. We 
determine the multifractal scaling properties through wavelet 
based detrended fluctuation analysis [17]. 



Materials and Methods 

Set up: A loud speaker was connected to a function gener- 
ator, a cathode ray oscilloscope and a resistor. A piezoelectric 
film was set properly on the base of the speaker to prepare the 
platform for the ball to bounce [12, 18]. Input frequency of 
the function generator was set to 49Hz and was made to gen- 
erate a sinusoidal potential, V = Vosin(ut). Fig.l depicts the 
experimental set up and Fig. 2 shows a typical potential time 
series. The 'Denoised' version of the raw signal is shown be- 
low it. For each input voltage Vo, the frequency input uj is 
varied from 2bHz to 100Hz, in the step size of IbHz. This 
process was repeated for 10V, 20V, 30V and 40V input volt- 
ages. The potential fluctuations of the piezoelectric crystal 
due to the bouncing of the ball above it, is the time series to 
be analyzed in the following sections. 
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at high frequencies the system's A is found to be positive. 
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FIG. 1 : Schematic illustration of the experimental setup. 
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FIG. 2: The potential time series of the bouncing ball with the cor- 
responding filtered one shown below it. The signal is for forcing 
parameter values of 40V and 70Hz. 

The bouncing ball shows periodic behavior for smaller val- 
ues of the forcing parameters. For higher values of the poten- 
tial, the time series exhibits intermittency [18]. It shows short- 
time periodic behavior and time varying transient dynamics. 
As the forcing frequency increases with constant forcing volt- 
age, the ball shows non-linear behavior [14]. 



Variation of Lyapunov exponent -with 
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FIG. 3: Variation in the Lyapunov exponent with constant forcing 
frequency, is shown as a function of forcing potential. The horizontal 
line corresponds to zero value of the Lyapunov exponent. These A 
values are obtained for 70Hz, forcing parameter. 




FIG. 4: Variation in the Lyapunov exponent with forcing potential 
constant, is shown as a function of input frequencies. The horizontal 
line corresponds to zero value of the Lyapunov exponent. These A 
values are obtained for 20V, forcing parameter. 

Figs. 3 and 4 depict the Lyapunov exponents obtained for 
varying both the forcing parameters. We now explore in more 
detail the dynamics of the potential time series. This is done 
using Fourier method and Wavelet transformation. 



ANALYSIS OF THE BOUNCING BALL MOTION 
Lyapunov Exponent 

We start with the estimation of the Lyapunov exponent (A). 
For certain low values of forcing parameters, it is found to be 
negative. As is well known, negative Lyapunov exponents are 
characteristic of dissipative systems. These systems exhibit 
asymptotic stability. Hence for negative A, the orbit attracts to 
a stable fixed point or stable periodic orbit. The more negative 
the exponent, the greater the stability. We also find that the 
bouncing ball system shows periodic behavior for small forc- 
ing parameters [19]. For higher values of the forcing param- 
eters, A > is obtained. In this domain, the system exhibits 
chaotic dynamics. We note that, even for low potential input, 



Frequency Domain Fourier Analysis 

FFT spectrum analyzer in vibrating dynamics is well 
known to measure frequency response of vibration or shock 
data. Fourier transform works by decomposing a signal to 
complex exponential functions of different frequencies [20] : 



/oo 
x(t)e- 27lift dt 
-OO 



and 



/oo 
X(/)e 2 "^/, 
-oo 

where / and t are frequency and time variables. 



(1) 



(2) 



Window-based finite impulse response filter and zero-phase 
digital filter have been used to remove noise, corresponding to 
the high frequency fluctuations. Fig. 5 shows the FFT corre- 
sponding to the filtered signal, with its unfiltered counter part. 



FFT of the Filtered signal 



FFT of the Raw signal 
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FIG. 5: FFT of filtered and raw signals are shown in 5 a and 5b re- 
spectively. 



Phase-space Plot 




FIG. 6: Phase-space plot corresponding to the denoised version of 
the signal, is clearly showing periodic motion. 

It is noticed that the noise strength in this time varying 
physical system, increases with increasing magnitude of the 
forcing parameters. The higher the magnitude of the frequen- 
cies present in the signal, the longer time it persists. As men- 
tioned earlier, we have applied a filter on the signal, which re- 
moves some initial as well as final transients to clearly observe 
the periodic behavior in the phase-space [5]. Fig. 6 depicts the 
phase-space plot of the filtered signal, which clearly shows 
the periodic dynamics. The number of converging loops show 
that the system has multiperiodicity. 
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FIG. 7: The power of the maximum output frequencies are shown for 
each input voltage to the system. 
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FIG. 8: The power of the maximum output frequencies are depicted 
for each input frequency. 

We now study the output from the system in order to iden- 
tify the parameter values corresponding to the non-linear mo- 
tion. Figs. 7 and 8 represent the response of the system to vary- 
ing forcing parameters, with the other one kept constant. In 
all signals, except for the input parameter of 10V, the highest 
frequency is same as the respective forcing frequencies indi- 
cating linear behavior. For 10V forcing potential, the highest 
response frequency, output is found to be 75Hz and 50Hz for 
25Hz and 100Hz input respectively. For all other forcing pa- 
rameters, the output response behaves linearly. For high val- 
ues of forcing parameters, chaotic behavior is observed from 
the frequency response. The frequency response, gets mini- 
mized due to dissipation of the system and also gets amplified 
for smaller values of forcing parameters. This clearly shows 
non-linear behavior [14]. 



CUMSUM Plot 





FIG. 9: [9(a)] Cumulative sum of the normalized profile and [9(b)] 
log-log plot of the FFT the corresponding profile are shown for 40V 
and 55Hz, forcing parameters. 



We now analyze the self- similar behavior of the time se- 
ries. For this purpose, the 'random walk' is calculated from 
the cumulative sum of the normalized profile of the system 
[20]. Fig. 9a shows the resulting time series, that reveals the 
fluctuation characteristic including its periodic components. 
Fig. 9b depicts log-log plot of the FFT of the cumulative sum 
of the normalized profile, which exhibits power law behavior 
[21], with varying exponents according to the forcing parame- 
ters. The Hurst parameter is estimated from the Fourier power 
law analysis, based on the monofractal hypothesis and gives 
the power law exponent /3 which is related to the Hurst expo- 
nent H by |/5| = 2i7 + 1 [20, 22-23]. It is observed that there 
are more than one exponent in the power law domain. They 
are calculated as /?i = 2.5526 and /?2 ~ 1, implying that 
H « in the second domain for 40V and 5 5 Hz forcing pa- 
rameters. In the first domain, persistence is observed; whereas 
the second domain, corresponding to the high frequency range 
of 330Hz to 25000Hz, indicates the possible presence of A 
behavior [25-26]. This structure shows the presence of mul- 
tifractality in the system, which can be estimated using the 
method like MFDFA [17]. 

As is well known, fractals can be subdivided in parts, which 
are a reduced copies of the whole system. They are character- 
ized by their fractal dimension D [27], which is a non-integer 
number that quantifies the density of fractals in the metric 
space. The number D measures the degree of fractal boundary 
fragmentation i.e., irregularity over multiple scales. 

The fractal behavior of the signal can be inferred as follows: 
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(3) 



Wavelet Based Multifractal Detrended Fluctuation Analysis 
[WBMFDFA] 

Wavelet analysis employs one series of different wavelet 
functions, which may consist of step functions instead of 
sine and cosine functions, derivative of Gaussian functions 
or a multipeak Morlet wavelet, composed of a sine wave 
multiplied with a Gaussian window. Discrete wavelet basis 
set is composed of a scaling function, mother and daughter 
wavelets, and has a strictly finite extent [28]. 

A convolution product of the data sequence (a function 
f(x), where x is usually a time or space variable) with the 
scaled and translated version of the wavelet function as, ip(x), 
is known as Wavelet transform. Two parameters are used 
to perform scaling and translation; the scale parameter s 
stretches (or compresses) the wavelet to the required reso- 
lution, while the translation parameter b shifts the reviewing 
wavelet to the desired location [29] : 



(Wf)(s,b) = - (°° f(x)r(—)da 

S J-oo s 



(4) 



where s, b are real as s > 0. The wavelet transform acts 
as a microscope: it uncovers finer details while operating at 
smaller scales [16]. 



Discrete Wavelet Transform (DWT) 

We have used discrete wavelets belonging to the 
Daubechies family to remove local polynomial trends from 
the profile. Considering Xi(i = 1, 2, 3, ...., N) to be the time 
series of data length N, the 'profile' Y(t) is the cumulative 
sum of series after subtracting the mean: 



Y(i) 



£ 
t = i 



x t — < X >\ 



1,2,3, N (5) 



Next, wavelet transform on the profile Y(t) has been im- 
plemented to separate the fluctuations from the trend. The 
fluctuation extraction involves a wavelet decomposition using 
the Db4 wavelet. For Db4, in addition to the defining relation 
f ipjidt = 0; the wavelets satisfy the vanishing moment con- 
dition f tipjidt = 0. This implies that for a linear signal of the 
type y = at + b, the wavelet coefficients will be zero. There- 
fore, this linear trend will be captured by the scaling function 
[16,21]. 

The profile series can be decomposed for the Db4 as 



In bouncing ball system, the value of D varies in between 
1 to 2, depending upon the forcing parameters. We now make 
use of wavelet transformation for a better understanding of the 
local dynamics. 



Y(t) 
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Y^ dijil>ij(t) 



(6) 



where, ipij (t) are wavelet functions of Db4 and fa (t) is the 
father wavelet. The coefficients q, d^ are known as the low 



pass and high pass coefficients. The profile is reconstructed 
at a particular level z, by taking only the low pass coefficients 
to extract the trend at level i. This trend is subtracted from 
the profile to obtain the fluctuations at each scale. Due to the 
convolution errors, these obtained fluctuations can have edge 
artifacts, which are removed by performing this fluctuation 
extraction on the reversed profile and taking the average. Af- 
ter that we subdivide the fluctuations in to n s non-overlapping 
segments, such that n s = [N/s], where the segment length s 
is related to the wavelet scale i by the number of filter coeffi- 
cients used in the wavelet. The q th order fluctuation function, 
F q (s) is obtained as follows [17]: 



Double-logarithmic plot of the Fluctuation function 



F q (s) 
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(7) 



Here q can take both positive and negative integral values. 
If the time series possesses fractal behavior, then F q (s) ex- 
hibits a power-law scaling, 



F q (s) 



Mq) 



(8) 



For q < 0, h(q), the generalized Hurst exponent, captures 
the scaling properties of the small fluctuations, whereas for 
q > 0, it captures the large fluctuations, q = 0, leads to 
divergence of the scaling exponent. In that case, logarithmic 
averaging has to be used to find the fluctuation function; 




FIG. 10: Log-log plot for the fluctuation function F q (s) vs. s for 
various values of q. This plot is for a signal having 40V and 70Hz as 
forcing parameters. 



Fig. 10 depicts the Double -logarithmic plot of the fluctua- 
tion function F q (s), for q between -10 to 10 in the ascending 
order. 

Local Hurst exponent variation with the order of moments 



q (order of moments) 
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For mono-fractal time series, the h(q) values are indepen- 
dent of q. For multifractal time series, h(q) depends on q. The 
correlation behavior is characterized from the Hurst exponent 
(H = h(q = 2)), which varies from < H < 1 [23-24]. 



RESULT 



Multifractal behavior of the system 



FIG. 1 1 : The multifractal behavior of the spectrum of the bouncing 
ball is shown. The dependence of the generalized Hurst exponent 
h(q), as a function of the order of the moment (q), is shown. The 
signal is for 70Hz and 40V as forcing parameters. 



Table. 1 shows the Hurst exponents obtained through DWT 
for different forcing parameters. 



TABLE I: Hurst exponents obtained for potential time series using 
Db4 wavelet 



Hurst exponents for varying forcing parameters 


Frequency (Hz) 


Potential (V) 


10 20 30 40 


25 


0.1981 0.0543 0.1090 0.5153 


40 


0.4652 0.2609 0.5429 0.6644 


55 


0.5806 0.5999 0.6778 0.7083 


70 


0.6286 0.6613 0.7116 0.7338 


85 


0.2860 0.7235 0.7173 0.7334 


100 


0.3466 0.6888 0.7171 0.7344 



The tabulated data explicitly shows that, the Hurst expo- 
nents depend on the forcing parameters. Some values of the 
Hurst exponent decrease for high forcing parameter with the 
other one at lower values. It increases with the high values 
of both the forcing parameters and vice versa. None of these 
multifractal trends could be gleaned using the Fourier analy- 



sis, where one can only qualitatively predict the overall nature 
of the self-similar behavior. Figs. 10 and 11 depict the com- 
mon nature of the Fluctuation functions and the generalized 
Hurst exponent respectively. It is found that for high values of 
the forcing parameters, our system shows persistence. 

We now proceed to analyze the periodic behavior through 
continuous wavelet transform (CWT). Both non- stationary 
periodic dynamics and the stationary ones clearly manifested 
in the scalogram of the CWT. We have used Morlet wavelet 
having optimal time-frequency localization. For a clearer 
identification of the periodic components, we removed the 
high frequency fluctuations through the Db4 basis. We have 
checked that the periodic components are not influenced in the 
process [26-27]. 



29]. 



3D plot showing wavelet coefficients from scales 1-2500 
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FIG. 12: 3D plot of CWT coefficients with time and scale axis, with 
40V and 70Hz as forcing parameters. Periodic behavior can be ob- 
served from these coefficients at certain scales. 



Semi-log plot of wavelet power 




FIG. 14: Semilog plot of wavelet power, summed over all time at dif- 
ferent scales. It can be inferred that, the ball has periodicity of 256ms 
and 579ms, corresponding to 40V and 70Hz forcing parameters. 



a)Bouncing ball motion 




FIG. 15: Wavelet power spectrum with the cone of influence (15b) 
for the signal (15a) having 40V and 70Hz as forcing parameters. The 
wavelet analysis can clearly distinguish the two dominant periods of 
the time series with the confidence level of 95%. Image colors are a 
representation of the Wavelet Power Spectrum (WPS) normalized by 
variance. The x-axis shows time; the logarithmic y-axis shows the 
wavelet scale. 
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From Figs. 13 and 14, it is evident that there are two dom- 
inant periods. Fig. 15 depicts the wavelet power spectrum in- 
cluding the global wavelet spectrum and the average variance 
with time. The inset in Fig. 15b clearly shows the presence 
of two periods in the potential time series. In the regions of 
higher power (the white regions), the dark spots inside the 
cone of influence, correspond to the periods of 256ms and 
579ms with 95% confidence level [29-30]. These periods 
compare well with the ones obtained through FFT earlier. 



CONCLUSION 



FIG. 13: Scalogram showing energy for each wavelet coefficients, 
with 40V and 70Hz as forcing parameters. 

Morlet wavelet clearly brings out different periodic compo- 
nents in the data sets, as seen in the scalograms given in Fig. 12 
in a 3D plot. Fig- 13 shows the scalogram of the potential time 
series, presenting the energy of the wavelet coefficients [28- 



In conclusion, we have studied the dynamics of a bounc- 
ing ball, a simple driven system exhibiting complex dynamics, 
which is controlled by the potential and frequency of the driv- 
ing source. Depending on the driving parameters, it shows 
self- similar and periodic behavior, which are studied using 
the Fourier domain analysis and Wavelet transform. It is ob- 
served that the potential time series shows a long term switch- 



ing between high and low frequency values in adjacent pairs. 
The characterization of the self- similarity revealed multifrac- 
tal behavior. Wavelets belonging to the Daubechies family 
optimally removed the local trends of the time series, from 
which the multifractal character has been reliably extracted. 
The continuous Morlet wavelet localized the periodic compo- 
nents, which revealed the linear and non-linear behavior. It is 
quite evident that the bouncing ball system provides an excel- 
lent experimental tool to study nonlinear dynamics and chaos. 
The apparatus can be used to study the complex chaotic or pe- 
riodic motions, depending on the parametric conditions of the 
ball. 

In future, we would like to analyze the i dynamics in more 
detail. We would also like to model the observed dynamics 
and explore the possibility of using this system to generate 
different random matrix ensembles [31-35]. 
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